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Abstract. By means of two- and three-dimensional particle-in-cell simulations, we 
investigate the process of driven magnetic reconnection at the termination shock 
of relativistic striped flows. In pulsar winds and in magnctar-powered relativistic 
jets, the flow consists of stripes of alternating magnetic field polarity, separated by 
current sheets of hot plasma. At the wind termination shock, the flow compresses 
and the alternating fields annihilate by driven magnetic reconnection. Irrespective 
of the stripe wavelength A or the wind magnetization a (in the regime <r > 1 of 
magnetically-dominated flows), shock-driven reconnection transfers all the magnetic 
energy of alternating fields to the particles, whose average Lorentz factor increases by 
a factor of a with respect to the pre-shock value. In the limit A/ (r^cr) >• 1, where is 
the relativistic Larmor radius in the wind, the post-shock particle spectrum approaches 
a flat power-law tail with slope around —1.5, populated by particles accelerated by the 
reconnection electric field. The presence of a current-aligned "guide" magnetic field 
suppresses the acceleration of particles only when the guide field is stronger than the 
alternating component. Our findings place important constraints on the models of 
non-thermal radiation from Pulsar Wind Nebulae and relativistic jets. 
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1. Introduction 

Pulsar Wind Nebulae (PWNe) are bubbles of synchrotron-emitting plasma powered by 
the relativistic wind of rotation-powered pulsars. Their broadband spectrum is produced 
by a nonthermal population of electron-positron pairs (hereafter, simply "electrons"), 
presumably accelerated at the wind termination shock, where the momentum flux of 
the relativistic pulsar wind is balanced by the confining pressure of the nebula. The 
observed radio spectrum (F u oc z/~ - 25 for the Crab Nebula; Bietenholz et al. 1997) 
implies a distribution of emitting electrons of the form dN/dE oc E~ p , with a flat slope 
p ~ 1.5. To explain the radio through optical emission of the Crab Nebula, the power 
law of shock-accelerated electrons should span at least three decades in particle energy 
(Lyubarsky 2003). Flat electron spectra with slopes p < 2 below GeV energies are also 
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required to model the radio emission of hotspots in radio galaxies (Stawarz et al. 2007) 
and X-ray observations of luminous blazar sources (Sikora et al. 2009). 

A flat power-law spectrum with index 1 < p < 2 is not expected from the standard 
theory of Fermi acceleration in relativistic shocks, which normally yields slopes p > 2 
(Achterberg et al. 2001, Keshet & Waxman 2005, Sironi & Spitkovsky 2009, Sironi & 
Spitkovsky 2011&). An acceleration mechanism capable of generating flat spectra in 
PWNe was discussed by Lyubarsky (2003) and studied by Sironi & Spitkovsky (2011a), 
under the assumption that the flow upstream of the termination shock consists of 
alternating stripes of opposite magnetic polarity, separated by current sheets of hot 
plasma (from now on, a "striped wind"). For obliquely-rotating pulsars, this is the 
configuration expected around the equatorial plane of the wind, where the sign of the 
toroidal field alternates with the pulsar period. A similar geometry is invoked in the 
proto-magnetar model of gamma-ray bursts, if the striped structure of the equatorial 
wind is preserved when the flow gets redirected along the polar jet (Metzger et al. 2011). 
When passing through a shock (e.g., the termination shock of pulsar winds), the flow 
compresses and the alternating fields annihilate by shock-driven reconnection, which 
may produce the flat electron spectrum required by the observations. 

The physics of magnetic reconnection can be captured self-consistently only 
by means of multi-dimensional particle-in-cell (PIC) simulations. Fully-kinetic PIC 
simulations provide a powerful tool to explore the microphysics of collisionless plasmas, 
since they can capture from first principles the fundamental interplay between charged 
particles and electromagnetic fields. In the context of relativistic magnetic reconnection 
in pair plasmas, most studies have explored the process of undriven reconnection 
(Zenitani & Hoshino 2001, Jaroschek et al. 2004), where field annihilation is initiated by 
a transient seed perturbation to an otherwise stable current sheet. As discussed above, 
this is not the setup expected at the termination shock of relativistic striped winds. 
There, it is the shock-compression of the flow that steadily drives regions of opposite 
magnetic field polarity toward each other, causing reconnection. Recent experimental 
and numerical studies of non-relativistic plasmas have shown that driven reconnection 
is much faster than the undriven process (Fox et al. 2011). 

In this work, we explore via multi-dimensional fully-kinetic PIC simulations the 
physics of driven reconnection at the termination shock of a striped relativistic electron- 
positron wind. We find that the alternating fields are completely dissipated upon 
compression by the shock, and their energy is transferred to the particles, regardless 
of the properties of the flow. Broad particle spectra with slopes 1 < p < 2 are a 
common by-product of shock-driven reconnection, but the extent of the power-law tail 
depends on the wind magnetization and the stripe wavelength. 

2. Simulation Setup 

We use the 3D electromagnetic PIC code TRISTAN-MP (Buneman 1993, Spitkovsky 
2005) to study the termination shock of a relativistic striped wind. PIC codes can model 
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astrophysical plasmas in the most fundamental way, as a collection of charged macro- 
particles that are moved by integrating the Lorentz force. The currents deposited by the 
macro-particles on the computational grid are then used to advance the electromagnetic 
fields via Maxwell's equations. The loop is closed self-consistently by extrapolating the 
fields to the macro-particle locations, where the Lorentz force is computed. 

We mainly utilize 2D computational domains in the xy plane, with periodic 
boundary conditions in the y direction. In §6 we compare 2D and 3D simulations 
and show that 2D runs can capture most of the relevant physics. For both 2D and 
3D domains, all three components of particle velocities and electromagnetic fields are 
tracked. The shock is set up by reflecting a magnetized electron-positron flow from a 
conducting wall located at x = (Fig. 1). The interaction between the incoming beam 
(that propagates along — x) and the reflected beam triggers the formation of a shock, 
which moves away from the wall along +x (Fig. 1). The simulation is performed in the 
"wall" frame, where the post-shock plasma is at rest. The incoming stream is injected 
with bulk Lorentz factor 70 = 15, but we have explored a wide range of 70, finding the 
same results, modulo an overall shift in the energy scale (Sironi & Spitkovsky 2011a). 

The upstream flow carries a strong magnetic field of intensity B , oriented along the 
y direction. The choice to initialize the field in the simulation plane (i.e., oriented along 
y, as opposed to z) is motivated by the agreement between our 2D simulations with in- 
plane fields and 3D experiments, as we show in §6. The field alternates with wavelength 
A. We vary A from A = 20c/w p up to A = 1280 c/cjp, where c/co p = a/ 'Jomc 2 / 4ire 2 n C Q 
is the relativistic plasma skin depth in the wind. Here, n c o is the density of particles 
in the wind. The thickness of each current sheet is a few plasma skin depths (wide 
enough so that the process of undriven reconnection is suppressed in our simulations), 
but our results are not sensitive to this parameter, provided that the contribution of 
current sheets to the particle and energy flux in the wind is negligible (see Appendix 
A). We parameterize the field strength B via the so-called magnetization parameter 
a = Bq 1 1 An^omn c oc 2 . We vary a from 10 up to 100. Although the magnetic field strength 
in the wind is always B , the wavelength-averaged field (B y )\ can vary from zero up 
to Bq, depending on the relative widths of the regions of positive and negative field 
(lower panel in Fig. 1). In pulsar winds, one expects {B y )\ = only in the equatorial 
plane (where the stripes are symmetric), whereas \(B y )\\/ Bq — > 1 at high latitudes 
(upper panel in Fig. 1). We choose a = 2(B y )\/(B + \ (B y ) x \) = 0.1 (corresponding to 
(B y )x/B ~ 0.05), to represent the behavior of the shock in the vicinity of the equatorial 
plane. Sironi & Spitkovsky (2011a) have explored the full range < a < 1. We remark 
that a net field {B y )\ 7^ does not play the role of a guide field, which in our geometry 
would correspond to a uniform component of magnetic field in the z direction (aligned 
with the particle velocity in the current sheets). The guide field is expected to be 
negligible in pulsar winds, but it may be dominant in gamma-ray burst and blazar jets. 
We explore the effects of a guide field on the physics of shock-driven reconnection in §7. 

We choose the relativistic electron (or positron) skin depth in the cold wind c/u p 
such that the smallest scale in the system, which for a > 1 is the relativistic Larmor 
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Figure 1. Upper panel: poloidal structure of the striped pulsar wind, according 
to the solution by Bogovalov (1999). The arrows denote the pulsar rotational axis 
(along fi, vertical) and magnetic axis (along /x, inclined). Within the equatorial 
wedge bounded by the dashed lines, the wind consists of toroidal stripes of alternating 
polarity, separated by current sheets (dotted lines) . Lower panel: simulation geometry. 
For 2D runs, the simulation domain is in the xy plane, oriented as shown in the upper 
panel. The incoming flow propagates along — x, and the shock moves away from the 
reflecting wall (located at x = 0) toward +x. The magnetic field lies in the simulation 
plane along the y direction, and its polarity alternates with wavelength A. A net 
stripe-averaged field (B y )\ > is set up by choosing red stripes (where B y = +B ) 
wider than the blue stripes (where B y = —Bo). For the pulsar wind sketched above, 
(B y )\ > is realized below the equatorial plane. 

radius r L = (c/u> p )/ y/a, is resolved with at least a few computational cells. | We usually 
choose rL = 3 cells, but we have verified that our results do not change when doubling 
the resolution. Our computational domain is typically ~ 400c/u p wide (along the y 
direction). As we show in Appendix B, a large box is essential for the consistency of our 
findings. Each computational cell is initialized with two electrons and two positrons, 
but we have performed limited experiments with a larger number of particles per cell 
(up to 32 per species), obtaining essentially the same results. 

3. Shock Structure 

In this section, we investigate the steady-state structure of a relativistic shock 
propagating in a striped electron-positron wind. We adopt 7 = 15, a = 10, and 

X The Debye length of the cold upstream wind is also resolved with a few computational cells, to avoid 
the finite grid instability and artificial heating. 
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Figure 2. Internal structure of the flow at uj p t — 3750, for A = 640c/u;p, 
a = 10, and a = 0.1. See the supplementary material for an associate movie 
(ld_shock_structure.mov) following the complete evolution of the shock. At w p i = 3750, 
the hydrodynamic shock is located at x ~ 1000 c/to p , and the MHD shock at 
x ~ 2600c/w p , as indicated by the vertical dotted blue line. As a function of the 
longitudinal coordinate x, we plot: (a) y-averaged particle number density n, in units 
of the upstream value; (b) y-averaged magnetic energy fraction cb = B 2 / 8irjomn c oC 2 
(black line), and y-averaged magnetic energy per particle es/n (red line); (c) y- 
averaged magnetic field B y (black line) and electric field E z (red line), normalized 
to the upstream magnetic field Bq; (d) mean kinetic energy per particle, in units of 
the bulk energy at injection; (e) (f) x — 7/3 x and x — 7/3 z positron phase space. 



(B y )x/B ~ 0.05 (corresponding to a = 0.1) as our fiducial values for the upstream 
bulk Lorentz factor, magnetization, and strength of the stripe-averaged field. To clarify 
the physics of shock-driven reconnection, we employ a relatively large value for the 
stripe wavelength (A = 320c/co> p or A = 640c/co> p ), and we neglect the guide field. In 
§5 we discuss the dependence of our results on the stripe wavelength and the wind 
magnetization, and in §7 we comment on the effects of a nonzero guide field. 

The longitudinal profile of the shock transition region, along the direction x of 
shock propagation, is presented in Fig. 2, whereas Fig. 3 shows the flow properties in 
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the xy plane of the simulation, zooming in on a region around the shock (as delimited by 
the vertical dashed red lines in Fig. 2(a)). The steady-state structure of the flow shows 
the presence of two shocks. The main shock (which we shall call the "hydrodynamic 
shock" from now on, for reasons that will become clear below) corresponds to the jump 
in density occurring at x ~ 1000 c/u p in Fig. 2(a). At x ~ 2600 c/u p , well ahead of the 
hydrodynamic shock, the incoming flow crosses a fast magnetohydrodynamic (MHD) 
shock, whose location is highlighted in Fig. 2 by a vertical dotted blue line. 

At the fast shock, the wind is decelerated and compressed. The compression 
induced by the passage of the fast shock through a given current sheet triggers the 
onset of magnetic reconnection. Small-scale islands develop inside the current sheet 
as a result of the tearing-mode instability. The growth and evolution of reconnection 
islands downstream from the fast shock (x < 2600 c/u p ) is clearly shown in the 2D plots 
of density and magnetic energy of Fig. 3(a) and (b), respectively. A magnetic X-point 
exists in between each pair of neighboring islands, where field lines of opposite polarity 
break and reconnect, and the magnetic energy is transferred to the particles. Particles 
flow into the X-points and are accelerated by the reconnection electric field (see the 
green and blue regions in the 2D plot of e B — e E of Fig. 3(c)), while being advected by 
the reconnected magnetic field into the closest island. Magnetic islands act as reservoirs 
of particles (Fig. 3(a)) and particle energy (the high-energy regions of Fig. 3(d), colored 
in red, match very well the overdense islands of Fig. 3(a)). As the flow recedes from the 
fast shock, the small-scale islands created inside each current sheet grow and coalesce 
to form bigger islands. This temporal evolution corresponds, in the snapshot of Fig. 3, 
to a pattern of bigger (and so, fewer) islands for current sheets that are farther behind 
the fast shock. Reconnection proceeds with distance downstream from the fast shock, 
up to the point where reconnection islands fill the entire region between neighboring 
current sheets. Now, the striped structure of the flow is erased, and a hydrodynamic 
shock forms (located at x ~ 1000 c/uo v in Figs. 2 and 3). 

Behind the hydrodynamic shock, the flow comes to rest (in the wall frame of our 
simulations), and the particle distribution is isotropic in three dimensions (compare 
longitudinal and transverse phase spaces in Fig. 2(e) and (f), respectively). The post- 
shock number density approaches n& ~ 4n u (Fig. 2(a), where n u is the density ahead 
of the fast shock), and the shock velocity is f3 s \ 1 ~ 1/3, in agreement with the jump 
conditions of a relativistic unmagnetized 3D plasma. In fact, most of the energy per 
particle downstream from the hydrodynamic shock is in kinetic form (as opposed to the 
dominant electromagnetic component of the incoming striped flow). This becomes clear 
when comparing the profile of the mean kinetic energy per particle (Fig. 2(d)) with the 
ID plot of the mean magnetic energy per particle (red line in Fig. 2(b)). f Behind the 
hydrodynamic shock, the average kinetic energy per particle reaches (7)/7o — cr+1 ~ 11, 

\ We point out that in the incoming flow electric and magnetic fields equally contribute to the energy 
balance, whereas the red line in Fig. 2(b) only includes magnetic fields. However, downstream from 
the hydrodynamic shock, the flow is at rest (in the simulation frame), and electric fields vanish (see 
the red line in Fig. 2(c)). 
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Figure 3. Internal structure of the flow at uj p t = 3750, for A = 640 c/w p , a = 10, 
and a — 0.1, zooming in on a region around the shock, as delimited by the vertical 
dashed red lines in Fig. 2(a). See the supplementary material for an associate movie 
(2d_shock_structure.mov) following the complete evolution of the shock. At u) p t — 3750, 
the hydrodynamic shock is located at x ~ 1000 c/uj p , and the fast MHD shock at 
x ~ 2600 c/w p . We show the 2D plots of: (a) particle number density, with contours 
showing the magnetic field lines; (b) magnetic energy fraction e^; (c) cb — £e, where 
6e = -E 2 /87T7o?7m c oC 2 ; (d) mean kinetic energy per particle. Panels (e) and (f) show 
the particle energy spectra, respectively outside (e) or inside (f ) of current sheets. The 
color of each spectrum matches the color of the corresponding arrow at the bottom of 
panel (d), showing where the spectrum is computed. The dotted line in panel (e) is 
a Maxwellian distribution with the same average energy as the downstream particles; 
the dashed lines in panels (e) and (f) indicate a power law with slope p = 1.4. 



as expected in the case of full dissipation of magnetic fields. Correspondingly, the mean 
magnetic energy per particle drops to zero (red line in Fig. 2(b)). Behind the shock, the 
only residual magnetic field comes from shock-compression of the stripe-averaged field 
(-B y ) a (this is best seen in experiments with a smaller stripe wavelength). 

Finally, we remark that, as apparent in Fig. 3, the tearing mode instability, which 
plays an essential role for the dynamics of the flow, can be captured correctly only with 
multi-dimensional simulations. Indeed, for the parameters employed in this section, 
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the one-dimensional model by Petri & Lyubarsky (2007) would predict negligible field 
dissipation, in sharp contrast with our findings. As we show in Appendix B, this clearly 
emphasizes the importance of multi-dimensional physics for our understanding of shock- 
driven magnetic reconnection. 

4. Particle Spectrum and Acceleration 

The particle energy spectrum at different locations through the flow, as marked by 
arrows at the bottom of Fig. 3(d), is shown in Fig. 3(e) and (f). Panel (e) shows the 
particle distribution in the cold wind, whereas panel (f) focuses on the hot plasma 
within current sheets. Downstream from the hydrodynamic shock, the striped structure 
is completely erased, and no difference persists between current sheets and cold wind. 
Here, the particle spectrum (black line in both Fig. 3(e) and (f)) is in the form of a 
flat power-law tail (with spectral index p ~ 1.4, dashed line) extending from 7 min ~ 30 
to 7 max ~ 500, where it cuts off exponentially. For comparison, a 3D Maxwellian with 
the same average Lorentz factor (7) = 70 (a + 1) ~ 165 is plotted in Fig. 3(e) as a 
dotted line, to show that the actual particle spectrum is much broader than a thermal 
distribution. By following the flow from the fast to the hydrodynamic shock, we now 
clarify how such a particle spectrum is generated. 

Ahead of the fast shock, the spectrum at x ~ 2850 c/ou p (yellow line in Fig. 3(f)) 
shows, in addition to the cold plasma in the striped wind, the hot particles that are 
initialized in the current sheet. Behind the fast shock (red line for x ~ 2200 c/u p , 
brown line for x ~ 1350c/cu p ), the particle spectrum consists of two components: a low- 
energy peak, reminiscent of the particle distribution far from current sheets (compare 
with the spectrum measured in the wind, dark blue line in Fig. 3(e)); and a high-energy 
component, that grows in number and extends to higher and higher energies, as the flow 
propagates toward the hydrodynamic shock. This component is populated by particles 
that were initially outside the current sheet, and have entered the sheet in the course 
of the reconnection process. By the time the flow enters the hydrodynamic shock, the 
current sheet occupies the entire flow, which explains the smooth transition from the 
spectrum within current sheets to the distribution downstream of the hydrodynamic 
shock (yellow through black lines in Fig. 3(f)). 

We now investigate in more detail the physics of particle acceleration within current 
sheets. In Fig. 4, we follow the trajectories of a representative sample of positrons, 
extracted from the simulation of a striped wind with A = 320 c/u p (but our conclusions 
hold for all wavelengths). All the selected positrons are initially in the cold wind (i.e., 
far from current sheets), at roughly the same x-location. They encounter the closest 
current sheet at u p t = 1008, as shown in Fig. 4(b) together with the corresponding 2D 
plot of the particle number density. Starting from this time, their histories can diverge 
significantly, depending on their y-location at the moment of interaction with the current 
sheet. The particles that will end up with relatively low energies (7 < 370, black lines in 
panel (a)) are found at u p t = 1008 exclusively around magnetic islands (their locations 
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Figure 4. Energy evolution of a sample of representative positrons interacting with 
a current sheet, in a striped flow with A = 320 c/oj p , a = 10, and a — 0.1. See 
the supplementary material for an associate movie (acceleration_general.mov) following 
the complete trajectories of the selected particles. Panel (a): Energy histories of the 
selected positrons, with lines of different color depending on their final energy (black if 
7 < 370, blue or red if 7 > 670) and their y-location at uj p t — 1008 (blue if y < 100 c/w p , 
red if y > 100c/w p ). Panel (b): Magnetic field lines (white contours) superimposed 
on the 2D plot of particle number density, at time uj p t = 1008 (marked as a vertical 
dashed line in panel (a)); depending on the final particle energy, the locations of the 
selected positrons are shown as yellow filled (if 7 < 370) or open (if 7 > 670) circles. 

are indicated with filled circles in panel (b)). At later times, they will get trapped on 
the outskirts of the growing islands, without appreciable changes in energy. 

On the contrary, the positrons that will eventually reach high energies (7 > 670, 
red and blue lines in Fig. 4(a)) are concentrated at u p t = 1008 in the vicinity of 
magnetic X-points (the particle locations are shown with open circles in Fig. 4(b)). 
Starting from u p t = 1008 (vertical dashed line in panel (a)), the energies of these 
positrons grow simultaneously and explosively (red and blue lines in Fig. 4(a)), as 
they are accelerated along z by the reconnection electric field parallel to the X-line 
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(Zenitani & Hoshino 2001, Jaroschek et al. 2004, Lyubarsky & Liverts 2008). Particles 
pre-accelerated at X-points continue to gain energy by compression (Lyubarsky & 
Liverts 2008, Giannios 2010), as they are advected along y from the X-point into the 
closest island. No major systematic changes of energy are observed after the particles 
get trapped within magnetic islands. Different X-points may provide different energy 
gains, as shown by the red (respectively, blue) curves in panel (a), which refer to high- 
energy positrons that get accelerated at the X-point above (respectively, below) the 
small magnetic island in the middle of the current sheet. In summary, the energy gain 
of a given particle is determined on the one hand by the strength of the reconnection 
electric field at the X-point, and on the other hand by the time available for acceleration, 
as the particle flows from the X-point into the closest island. As the flow propagates 
downstream from the fast shock, the typical distance between an X-point and the closest 
island increases, allowing particles to be accelerated to higher energies. This explains 
why the upper cutoff in the particle spectra of Fig. 3(f) shifts to higher energies for 
current sheets farther behind the fast shock. 

The validity of our conclusions may be questioned if the structure of X-points in a 
realistic 3D scenario is much different from the 2D picture analyzed here. In particular, 
in 3D we expect the current sheet to be folded, by effect of the so-called drift-kink 
instability (Zenitani & Hoshino 2005, Zenitani & Hoshino 2007). The characteristic 
wavelength of the kink mode in the xz plane will introduce a different length scale, 
that could potentially compete with the distance between X-points and islands (in 
the xy plane) in constraining the maximum energy of accelerated particles. In §6 we 
address this important point, and show that the latter constraint is more restrictive 
than the former, at least for the parameters explored in this work. It follows that our 
2D simulations with in-plane fields can capture very well the main 3D properties of the 
shock, and in particular the physics of particle acceleration. 

Although we find that the reconnection electric field at X-points is the main agent of 
particle acceleration, we have verified, by tracing the orbits of a large number of particles 
extracted directly from the simulation, that a variety of other mechanisms can play a role 
in particle energization. In Fig. 5, we show the energy evolution of three representative 
positrons (panel (a)), together with the 2D density structure at the time when a given 
positron is being rapidly accelerated (panel (b) corresponds to the vertical red dashed 
line in panel (a), panel (c) to orange, panel (d) to yellow). We find that the initial 
energy gain seen in the red and yellow curves of Fig. 5(a) (at u p t ~ 600 and cu p t ~ 1250, 
respectively) is due to the reconnection electric field that the particles experience upon 
their first encounter with an X-point, as described above. In contrast, the acceleration 
seen in the orange curve at u p t ~ 1100 is powered by the anti-reconnection electric field 
at the X-point located in between a pair of merging islands, as discussed by Jaroschek 
et al. (2004) and Oka et al. (2010). In fact, Fig. 5(c) shows that the orange trajectory 
at u p t ~ 1087 is confined in between two coalescing islands. In summary, reconnection 
and anti-reconnection at X-points govern the first stages of particle acceleration. 

The late-time energization seen in the red and yellow curves of Fig. 5(a) (at 
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Figure 5. Energy evolution of three representative high-energy positrons, in a striped 
flow with A = 320c/u;p, a — 10, and a = 0.1. See the supplementary material for an 
associate movie (acceleration_special.mov). Panel (a): Energy histories of the selected 
positrons. Panels (b)-(d): 2D plots of particle number density, at the time when the 
three positrons in panel (a) are experiencing fast acceleration, as indicated by the 
vertical dashed lines in panel (a). The positron trajectory is plotted with the same 
color coding as in panel (a). Along the positron track, the filled circle indicates the 
particle location at the time of the corresponding 2D density structure, and the arrow 
shows the velocity direction. 



u p t ~ 850 and u p t ~ 1800, respectively) is due to processes not primarily related to the 
reconnection electric field. The positron shown in Fig. 5(b) is accelerated by a kind of 
first-order Fermi process within magnetic islands, as discussed by Drake et al. (2006). 
The selected particle gains energy by reflecting from the ends of contracting islands 
following coalescence (its trajectory is orbiting around two merging islands in Fig. 5(b)). 
In contrast, the positron in Fig. 5(d) is accelerated at u p t ~ 1800 via the standard first- 
order Fermi acceleration at shocks (Blandford & Ostriker 1978, Bell 1978, Blandford & 
Eichler 1987), by bouncing off a magnetic island (located in the middle of Fig. 5(d)) 
that is advected into the shock with the upstream flow. 
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Figure 6. Downstream particle spectrum for different values of A, in a flow with 
a = 10 and a = 0.1. The dotted red line is a Maxwellian with the same average energy 
as the spectrum colored in red (which refers to A = 20 c/ui p ); the dashed black line 
is a power law with slope p — 1.5. In the subpanel, the black line shows the average 
Lorentz factor of downstream particles (7) as a function of A (axis on the left), whereas 
the red line presents the dependence of 7max/7min on A (axis on the right). Here, 7 m j n 
is the location where r ydN/d'f peaks (i.e., where most of the particles reside), whereas 
7max is the Lorentz factor where ^dN/d^f peaks (i.e., where most of the energy lies). 



5. Dependence on the Wind Wavelength and Magnetization 

In this section, we investigate how the physical conditions in the wind can affect the 
structure of the shock and the physics of magnetic reconnection, with particular focus 
on the particle spectrum downstream from the hydrodynamic shock. We explore the 
dependence of our results on the stripe wavelength and the wind magnetization, for the 
case with no guide fields. We show that complete dissipation of the alternating fields 
(and transfer of field energy to the particles) occurs in all cases, but the width of the 
downstream particle spectrum depends on the properties of the wind. 

Fig. 6 shows how the post-shock energy spectrum changes for different stripe 
wavelengths A, keeping fixed the magnetization a = 10 and the stripe-averaged field 
(B y )\/B ~ 0.05 (corresponding to a = 0.1). In all cases, the value of the post-shock 
mean particle Lorentz factor (7) measured in our simulations (black line in the subpanel 
of Fig. 6) is consistent with full dissipation of the alternating fields, yielding (7) ~ 7 cr. 
Even though the mean kinetic energy per particle does not appreciably vary with A, the 
shape of the spectrum does change, with a clear tendency for broader spectra at longer 
stripe wavelengths (as indicated also by the ratio 7max/7min between upper and lower 
spectral cutoffs, red line in the subpanel of Fig. 6). 
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Figure 7. Downstream particle spectrum for different values of A and a, but such that 
the ratio \/{tl<j) ~ 22.5 is kept fixed (as well as a = 0.1). Here, \/{ri,o) is roughly 
the wind wavelength in units of the post-shock plasma skin depth. The comparison 
is performed at the same time, in units of the relativistic plasma frequency of the 
post-shock flow (~ ujp/y/a). Spectra are shifted along the ir-axis by 50 /a to facilitate 
comparison with the case a = 50. The black line in the subpanel shows the average 
downstream Lorentz factor (7) as a function of a. 

This trend can be easily understood by considering the structure of the flow, 
just upstream of the hydrodynamic shock. Here, the size of reconnection islands is 
constrained by the distance between neighboring current sheets, which is proportional to 
A. Shorter wavelengths will then result in more numerous islands of smaller size, whereas 
fewer but bigger islands will be present for longer A. Since an X-point exists in between 
each pair of neighboring islands (belonging to the same current sheet), the number of 
X-points per unit length (along the current sheet) will be larger for smaller A. For short 
wavelengths, most of the incoming particles will likely pass in the vicinity of one of the 
numerous X-points, thus gaining energy from the reconnection electric field. All the 
particles will share similar energy evolutions, which results in a narrow Maxwellian-like 
distribution. On the other hand, for long wavelengths the energy evolution of different 
particles can be extremely diverse (as discussed in §4), depending on how far they pass 
from the closest X-point. Most of the particles stay far from X-points and remain cold, 
but the particles that interact with an X-point are accelerated to high energies by the 
reconnection electric field. This results in a broad energy spectrum. 

We find that the threshold between short and long wavelengths depends on the 
stripe wavelength and the wind magnetization via the combination A/ (V^cr), namely the 
stripe wavelength measured in units of the post-shock plasma skin depth. In Fig. 7 we 
confirm that the main properties of the shock, and in particular the high-energy end 
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of the downstream particle spectrum, are relatively insensitive to variations in A or a, 
provided that the ratio A/(r^o") is kept constant. When varying A/(r£<r), a Maxwellian- 
like spectrum is obtained for A/(r£(x) < a few tens, whereas in the limit A/(r^(T) ^> 1 
the spectrum approaches a broad power-law tail of index 1 < p < 2, extending from 
7min ^ 7o up to 7 max ~ 7 o- 1/(2 ~ p) . 

The minor bump emerging at high energies in the spectra of Figs. 6 and 7 is 
populated by particles that are accelerated via shock-drift acceleration (SDA; Chen & 
Armstrong 1975, Webb et al. 1983, Begelman & Kirk 1990) at the hydrodynamic shock. 
They gain energy from the stripe-averaged motional electric field (E z ) x ^ (B y }\, while 
gyrating around the shock, before being advected downstream by the stripe-averaged 
magnetic field. Given the limited number of acceleration cycles, the SDA component 
does not extend in time to higher energies. So, SDA is not a promising candidate to 
produce broad power-law tails. 

6. Three-Dimensional Simulations 

The results presented so far are based on 2D experiments with magnetic field initialized 
in the plane of the simulations. This is the best configuration to study the tearing-mode 
instability, whose wavevector lies in the plane of the field, but it artificially inhibits 
the growth of cross-plane instabilities, which may also be important for the structure 
of the shock. In this section, we study the 3D physics of shock-driven reconnection, 
with particular emphasis on particle acceleration. The simulation setup parallels very 
closely what we described in §2, with the alternating field oriented initially along the y 
direction and no guide field. We adopt our fiducial values for a — 10 and a = 0.1, and 
we investigate A = 80c/w p and A = 160c/w p . For each value of the stripe wavelength, 
the extent of the simulation domain along y and z is chosen such that our results are 
not artificially affected by the periodicity of our boundary conditions (see Appendix B). 

Fig. 8 shows the structure of the shock at u p t = 990, for a striped wind of wavelength 
A = 160c/wp. In agreement with the results presented in §3 for 2D simulations, the 
structure of the flow in Fig. 8 shows the presence of two shocks, a fast MHD shock 
located at x ~ 630c/w p and a hydrodynamic shock at x ~ 250c/w p . The passage of 
the fast shock through the incoming current sheets initiates magnetic field annihilation. 
The dissipation of alternating fields proceeds as the flow propagates from the fast to 
the hydrodynamic shock, and little magnetic energy remains downstream from the 
hydrodynamic shock (see Fig. 8(b) at x < 250c/w p ), in agreement with our 2D results. 

As the flow recedes behind the fast shock, the structure of the current sheets is 
affected by the growth of two competing modes. In the xy plane of the magnetic field, 
the tearing-mode instability breaks the current sheet into a sequence of high-density 
islands (see Fig. 8(a)), separated by X-points where the field lines tear and reconnect. 
In the xz plane orthogonal to the field, the folding of the current sheet seen in Fig. 8 is 
governed by the drift-kink instability, driven by the current of fast-drifting plasmas in a 
thin sheet (Daughton 1998, Zenitani k Hoshino 2005, Zenitani & Hoshino 2007). The 
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Figure 8. Internal structure of the flow at u> p t — 990, from the 3D simulation of 
a striped wind with A = 160c/o;p, a = 10, and a = 0.1. See the supplementary 
material for associate movies (3d_density.mov and 3d_magnetic_energy.mov) that follow 
the complete evolution of the shock. The magnetic field is initialized in the xy plane, 
as shown by the black arrows in panel (a). The hydrodynamic shock is located at 
x ~ 250c/w p , and the fast MHD shock at x ~ 630c/u;p. We show the 3D plots of: (a) 
particle number density, in units of the upstream value, with color scale stretched to 
enhance contrast; (b) magnetic energy fraction t B = i? 2 /87r7 TOn c0 c 2 , with color scale 
stretched to enhance contrast. 

2D geometry employed in the previous sections was chosen to select the tearing mode, 
and suppress the drift-kink mode. 

The nature of the instability that dominates the process of field annihilation should 
leave an imprint on the resulting particle distribution. Zenitani & Hoshino (2007) have 
shown with PIC simulations that field dissipation due to the drift-kink instability does 
not result in nonthermal particle acceleration, the plasma just being heated. This is 
because the field lines remain straight, so that all particles gain energy at the same rate. 
Nonthermal particles are produced only by the tearing mode (Zenitani & Hoshino 2001). 
In our 2D simulations with in-plane fields, i.e., the geometry required to capture the 
tearing mode, we have shown that particles are accelerated to nonthermal energies 
by the reconnection electric field. They initially move along the z direction parallel 
to the reconnection field, and then they drift along y from a given X-point into the 
closest magnetic island. In 3D, one could argue that the folding of the current sheet 
introduced by the drift-kink instability in the xz plane (not resolved in 2D) may deflect 
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Figure 9. In a striped wind with A = 80c/cj p , a = 10, and a = 0.1, comparison 
of the downstream particle spectrum for different magnetic field configurations: 3D 
simulation (red) and 2D simulations with either in-plane (blue) or out-of-plane (green) 
magnetic fields. The black line in the subpanel shows the average downstream Lorentz 
factor (7) for the different cases. 



the accelerating particles out of the current sheet, thus suppressing their energization. 

In Fig. 9, we address this important issue, by comparing for A = 80c/u;p the post- 
shock particle spectrum of a 3D simulation (red line) to the results of 2D experiments 
with in-plane or out-of-plane magnetic fields (blue and green lines, respectively). The 
excellent agreement between 3D and 2D in-plane results (red and blue line, respectively) 
suggests that the physics of particle acceleration by shock-driven reconnection is 
captured extremely well by our 2D simulations with in-plane fields. § In turn, this 
implies that, at least for A = 80c/u p , the maximum energy of accelerated particles is 
constrained by the distance between X-points and islands (in the xy plane), rather than 
by the current sheet folding in the xz plane. 

Finally, as we discuss in §7, the tension force of a guide field (along z, in our 
geometry) can easily stabilize the drift-kink mode, with little or no effect on the tearing 
mode (Zenitani & Hoshino 2008). In the presence of a guide field, the 3D physics of 
shock-driven reconnection should be described very accurately by 2D simulations with 
alternating fields lying in the simulation plane. 

§ The agreement between red and blue lines in the high-energy bump at 7 > 2000 also suggests that 
the physics of SDA is correctly described by our 2D simulations with in-plane fields. 
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Figure 10. 2D structure of the particle number density (in units of the upstream 
value), for alternating fields lying either in the simulation plane (left column) or 
orthogonal to the plane (right column). We vary the guide field strength B g /B as 
indicated on the left. The figure refers to a flow with A = 320c/w p , a = 10, and 
a = 0.1, at time u p t = 1500. 



7. Dependence on the Guide Field Strength 

The results presented so far are obtained in the absence of a guide field, i.e., the 
component of magnetic field aligned with the particle current in the sheets. In this 
section, we investigate the effects of a guide field component on the structure of the 
termination shock. A guide field is not usually expected in the context of pulsar winds 
(but see Petri & Kirk 2005), but it may be important in gamma-ray burst and blazar jets, 
where it could be even stronger than the alternating component. We assume that the 
alternating fields are oriented either along or perpendicular to the xy simulation plane 
(so that the guide field is along z or y, respectively). In addition to the alternating fields 
of strength B , we initialize in the pre-shock flow a uniform guide field of intensity B g 
(with 0.0 < Bg/ Bq < 1.0), both within and outside the current sheets. 

The role of the guide field in the case of undriven reconnection has been investigated 
by Zenitani & Hoshino (2008) by means of linear analysis and PIC simulations. They 
found that the tearing mode is relatively insensitive to the strength of the guide field, 
whereas the relativistic drift-kink instability is stabilized by the tension force of the 
guide field if B 9 /Bq > 0.3. In this section, we perform a similar analysis for shock- 
driven reconnection at the termination shock of striped winds, focusing on the process 
of particle acceleration. For alternating fields lying in the plane of the simulations, the 
drift-kink mode is not resolved, and we can isolate the effects of the guide field on the 
tearing mode. Instead, if the alternating fields are orthogonal to the computational 
domain, we can study whether the guide field can stabilize the drift-kink mode. 

Fig. 10 shows the structure of the flow in the two cases (left column for in-plane 
alternating fields, right column for out-of-plane fields), for different strengths of the guide 
field, ranging from B g /B = 0.0 up to 0.3. The wavelength and magnetization of the 
striped wind are A = 320 c/u p and o = 10 in all cases. The left column shows that, in the 
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Figure 11. Downstream particle spectrum at w p i = 3000 for different values of the 
guide field strength B g /B , in a flow with A = 320c/wp, a = 10, a = 0.1, and 
alternating fields lying within the simulation plane. In the subpanel, the black line 
shows the mean Lorentz factor (7)770 of post-shock particles as a function of B 9 /Bq. 



case of in-plane alternating fields (so that the guide field is along z), the development of 
shock-driven reconnection and the resulting structure of the flow are nearly insensitive to 
the strength of the guide field. For all the values of B g we investigate (B g / B from 0.0 to 
0.3, from top to bottom), a fast shock propagates into the striped wind, compressing the 
incoming current sheets and triggering the tearing mode, which results in the magnetic 
islands seen in the density structure of Fig. 10 (left column). When the size of magnetic 
islands becomes comparable to the distance between neighboring current sheets, the 
striped structure of the flow is erased, and the main shock forms (at x ~ 400 c/u p in the 
left column of Fig. 10), as described in §3. The alternating fields annihilate completely 
at the shock, whereas the guide field is compressed, in agreement with MHD jump 
conditions. For stronger guide fields, the magnetic pressure in the post-shock plasma is 
larger, and the shock propagates at a higher speed. This explains why the location of 
the main shock shifts toward +x as the guide field strength increases from 0.0 to 0.3 
(panel (a) through (d) in the left column). 

In summary, the overall flow structure in the case of in-plane alternating fields - 
the configuration required to capture the tearing mode instability - is not appreciably 
changed by the presence of the guide field, up to B g /B = 0.3. This is confirmed by 
Fig. 11, where we compare the post-shock particle energy spectra for values of B g /B 
ranging from 0.0 to 1.0. We find that the particle spectrum departs from the result 
obtained in the absence of a guide field (black line) only for B g /B > 0.3 (yellow line 
for B g /B = 0.3, red line for B g /B = 1.0). For strong guide fields the rate of particle 
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Figure 12. Downstream particle spectrum at w p i = 3000 for different values of the 
guide field strength B g /B , in a flow with A = 320c/wp, a = 10, a = 0.1, and 
alternating fields orthogonal to the simulation plane. In the subpanel, the black line 
shows the mean Lorentz factor (7)/7o of post-shock particles as a function of B 9 /Bq. 

inflow into the X-points is smaller, since the plasma should be moving across the guide 
field lines. This results in fewer particles being accelerated by the reconnection electric 
field. This explains why the high-energy component in the particle spectra of Fig. 11 is 
suppressed for B g /B > 0.3. In turn, this corresponds to the decrease of mean particle 
Lorentz factor shown in the subpanel of Fig. 11, with respect to the value (7) ~ 7o(o" + 1) 
expected in the case of complete field annihilation. || 

So far we have discussed the effects of a guide field component on the development 
of the tearing instability, the only mode expected for alternating fields lying in the 
simulation plane. The right column of Fig. 10 shows instead the effect of the guide field 
in the case of alternating fields oriented along z, i.e., perpendicular to the simulation 
plane. For this geometry, the guide field is oriented along y. The flow structure in 
the absence of a guide field (B g /B = 0.0 in panel (e) of Fig. 10) shows a fast MHD 
shock located at x ~ 1400 c/u p , which compresses the incoming flow and triggers the 
drift-kink mode, manifested in the folding of the current sheets downstream from the 
shock. The right column of Fig. 10 clearly shows that the drift-kink mode is stabilized 

|| The different normalization of the spectral peak for B 9 /Bq = 1.0 (red line in Fig. 11) results from 
the fact that the post-shock flow is more magnetized, for extreme values of the guide field. It follows 
that the shock speed is higher, so the post-shock density is lower. 

% The speed of the shock is much larger here than in the case of in-plane fields with B g /Bo — 0.0 of 
Fig. 10(a). This results primarily from the different adiabatic index imposed by the 2D geometry of 
our simulations. If the alternating fields are in the simulation plane, the particle gyro-motion around 
the field will involve the third dimension (i.e., orthogonal to the plane), so that the post-shock plasma 



Shock-driven reconnection in relativistic striped winds 



20 



once the guide field exceeds B g / Bq ~ 0.1, in agreement with the conclusions by Zenitani 
& Hoshino (2008). For B g /B = 0.3 in Fig. 10(h), the current sheets downstream from 
the shock do not show any sign of the kink mode. Due to the suppression of the drift- 
kink mode for B g /Bo > 0.1, most of the post-shock energy stays in the magnetic field, 
and the plasma remains cold. This is clearly shown in the particle energy spectra of 
Fig. 12, for different strengths of the guide field. As the value of B g /B increases, the 
spectrum passes from a flat distribution with mean Lorentz factor ~ 7o(c + 1) (black 
line for B g /B = 0.0), as expected for efficient field annihilation, to a mono-energetic 
distribution peaking at 7 ~ 70 (red line for B g /B = 1.0), as in the case of negligible 
field dissipation. The subpanel of Fig. 12 shows that the suppression of the drift-kink 
mode, and the resulting inhibition of particle heating, occurs for guide field strengths 
as small as B g /B ~ 0.03. 

8. Summary and Discussion 

We have explored by means of 2D and 3D PIC simulations the internal structure 
and acceleration properties of relativistic shocks that propagate in an electron-positron 
striped wind, i.e., a flow consisting of stripes of alternating field polarity separated by 
current sheets of hot plasma. 

We find that a fast MHD shock propagates into the pre-shock striped flow, 
compressing the incoming current sheets and initiating the process of driven magnetic 
reconnection, via the tearing-mode instability. Reconnection islands seeded by the 
passage of the fast shock grow and coalesce, while magnetic energy is dissipated at 
X-points located in between each pair of islands. When reconnection islands grow so 
big to occupy the entire region between neighboring current sheets, the striped structure 
of the flow is erased, and a hydrodynamic shock forms. By this point, the energy stored 
in the alternating fields has been entirely transferred to the particles via shock-driven 
magnetic reconnection, so the post-shock fluid behaves like an unmagnetized plasma. In 
our 2D and 3D simulations, we find that complete annihilation of the alternating fields 
occurs irrespective of the stripe wavelength A or the magnetization a. 

The main agent of particle energization is the reconnection electric field, as particles 
drift from a given X-point into the closest island. Whether all particles have comparable 
energy gains, or only a few particles are accelerated to high energies, and the majority 
stay cold, depends sensitively on the properties of the wind. In the absence of guide 
fields, we find that the shape of the post-shock spectrum depends primarily on the 
combination X/{tl(t), where is the relativistic Larmor radius in the striped wind. 
For small values of X/ir^a) (< a few tens), the spectrum resembles a Maxwellian 
distribution with mean energy (7) ~ 7 a, where 70 is the bulk Lorentz factor of the pre- 

will be isotropic in 3D. In contrast, for out-of-plane alternating fields, the motion of the plasma in the 
presence of the field is confined to the 2D plane of the simulations, resulting in a 2D adiabatic index. 
The fact that the tearing mode is more effective than the drift-kink mode in converting field energy to 
particle energy further contributes to the difference in shock speed between panel (a) and (e). 
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shock flow. In the limit of very large values of \/{tlO') (> a few hundreds), the spectrum 
approaches a broad power-law tail with slope p ~ 1.5, extending from 7 min ~ 7 up to 
7 max ~ 7 (7 1 ^ 2_p - ) . Our results are nearly insensitive to the strength of the guide field B g , 
in the regime B 9 /Bq < 0.3, where Bq is the magnitude of the alternating component. 
For strong guide fields (i.e., B g /B > 0.3) the rate of particle inflow into the X-points is 
smaller, since the plasma should be moving across the guide field lines. This results in 
fewer particles being accelerated by the reconnection electric field, and the post-shock 
spectrum shifts to lower energies, everything else being fixed. 

The particle energy spectra extracted from our simulations can be used directly to 
interpret the radiative signature of astrophysical nonthermal sources, and in particular 
of PWNe. The radio spectrum of the Crab Nebula, the prototype of PWNe, requires a 
population of nonthermal particles with a flat spectral slope (p ~ 1.5), extending from 
10 2 MeV up to 10 5 MeV. Our results suggest that, for p ~ 1.5, the particle spectrum 
will extend over three decades in energy if the wind magnetization at the termination 
shock is a > 30. In addition, broad particle spectra are generated only if \/{ri,a) > a 
few tens. At the termination shock of the pulsar wind (R = Rts) we have 

~ Attk^^- , (1) 



r L a R TS 

where Rlc = is the light cylinder radius (Q is the pulsar rotational frequency), 
and k is the so-called multiplicity in the wind (i.e., the ratio of the actual density to 
the Goldreich & Julian's (1969) density). For the Crab, Rts — 5 x 10 8 R lc (Hester 
et al. 2002), and most available models estimate k ~ 10 4 — 10 6 (Bucciantini et al. 2011). 
Based on our findings, the resulting value of A/(rL<r) < 0.01 would yield a Maxwellian- 
like spectrum, at odds with the wide flat spectrum required by observations. If radio- 
emitting electrons are produced in the equatorial plane by shock-driven reconnection, 
a revision of the existing estimates of k is required. In this respect, we point out that 
the values of k quoted above are averages over latitude, and one cannot exclude that 
particle injection into the pulsar wind is highly anisotropic, with multiplicity as large as 
10 8 along the equatorial plane. This is not in conflict with the observed lack of gradients 
in the radio spectral slope of the Crab (Bietenholz & Kronberg 1992), since strong fluid 
motions downstream from the termination shock would quickly fill the entire nebula 
with the long-lived radio-emitting electrons (Del Zanna et al. 2004). Alternatively, the 
radio part of the spectrum may be produced at higher latitudes, where the termination 
shock is closer to the pulsar, and the ratio in eq. (1) becomes larger. In summary, our 
findings place important constraints on the magnetospheric physics of pulsars, and they 
provide physically-grounded inputs for models of nonthermal emission in PWNe. 
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Appendix A. Dependence on the Current Sheet Thickness 

The results presented in the main body of the paper have been obtained assuming that 
the thickness of current sheets in the wind is comparable to the relativistic skin depth, 
i.e., A ~ c/up. This value is large enough so that the process of undriven reconnection 
does not significantly modify the flow structure ahead of the termination shock. On 
the other hand, since 5 = 7rA/A <C 1, the current sheets only occupy a small fraction 
of the striped flow, and most of the pre-shock particles and energy are contributed by 
the material in the cold wind. In Fig. Al, we study the dependence of the post-shock 
particle spectrum on the current sheet thickness, by increasing the fractional width S 
from 0.01 (as adopted so far for A = 320c/cj p ) up to 1.0. 

cr = 10 A=320 c/cj p a = 0A u p t=3000 
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Figure Al. Downstream particle spectrum at ui p t = 3000 for different values of the 
current sheet width A, parameterized in terms of the fractional thickness 6 = irA/X, in 
a flow with A = 320 c/ui p , a — 10, and a = 0.1. We vary S from 8 = 0.01 (corresponding 
to A ~ c/a)p) up to S — 1.0 (corresponding to A ~ A/V, i.e., comparable to the half- 
wavelength of the stripes). The solid line in the subpanel shows the average downstream 
Lorentz factor (7), i.e., the post-shock mean kinetic energy per particle. The dotted 
line shows the pre-shock mean kinetic+ electromagnetic energy per particle, and the 
dashed line is the pre-shock mean electromagnetic energy per particle. The different 
spectra are normalized such that to contain the same number of particles. 

As shown in Fig. Al, for S < 0.2 the shape of the post-shock spectrum is insensitive 
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to the sheet thickness. In this case, both the evolution of the flow across the shock, and 
the physics of particle acceleration by reconnection, proceed in the same way as described 
in §3 and §4. The only difference is in the size of magnetic islands just downstream from 
the fast MHD shock, since the wavelength of the most unstable tearing mode scales 
linearly with the sheet thickness (Zenitani & Hoshino 2007). 

As 5 increases beyond 0.2, the high-energy end of the post-shock spectrum gets 
progressively depopulated, and it shifts to lower energies. At the same time, the post- 
shock mean kinetic energy per particle (solid line in the subpanel) decreases." 1 " However, 
this is not the result of a smaller efficiency in converting electromagnetic energy to 
particle energy via reconnection. Rather, it comes from the fact that, with increasing 
S, the mean kinetic-h electromagnetic energy per particle in the pre-shock flow becomes 
smaller, everything else being fixed (dotted line in the subpanel). This is because in our 
setup the current sheets are overdense by a factor of four with respect to the cold wind, 
so that the mean particle energy within current sheets is smaller by the same amount. 
In the subpanel, the agreement between the solid line (post-shock mean kinetic energy 
per particle) and the dotted line (pre-shock mean kinetic-h electromagnetic energy per 
particle) suggests that complete field dissipation occurs regardless of 5. However, as 5 
increases, a smaller fraction of the pre-shock energy is in electromagnetic form (dashed 
line in the subpanel), which corresponds to a smaller effective magnetization of the 
upstream flow. It follows that the reconnection electric field will be weaker, and particles 
will not be accelerated to energies as high as in the case of thin current sheets. 

Appendix B. Dependence on the Transverse Size of the Simulation Box 

In Fig. Bl we investigate the dependence of our findings on the transverse size L y of 
the computational domain, for A = 640c/w p , a = 10, and a = 0.1. The agreement 
between the black curve [L y = 400c/u p ) and the red line [L y = 800c/w p ) suggests that 
our results become insensitive to the transverse size of the box (i.e., they converge with 
respect to L y ), for boxes larger than 400c/w p . Generally speaking, we find that a good 
criterion for convergence is L y > A/2.* When this condition is fulfilled, the growth of 
magnetic islands behind the fast shock is not artificially inhibited by the periodicity of 
our boundaries in the y direction. Magnetic reconnection then proceeds up to the point 
when islands from two neighboring current sheets will merge, which happens when their 
size is ~ A/2 (for a < 0.1). 

On the other hand, for L y < 300c/o;p the growth and coalescence of reconnection 
islands within a given current sheet artificially stops when only one island is left in 

+ Strictly speaking, the solid line also includes the contribution of the rest-mass energy. In this section, 
we simply call "kinetic energy" (as opposed to "electromagnetic energy") the sum of the kinetic and 
rest-mass energy. 

* We point out that this criterion is much more constraining than the requirement that the fastest 
growing tearing mode can fit in the width of the computational box. For relativistic reconnection, 
the most unstable wavelength is comparable to the current sheet thickness (Zenitani & Hoshino 2007), 
which in our setup is a few skin depths, much smaller than the stripe wavelength A = 640c/cj p . 
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Figure Bl. Downstream particle spectrum at uj p t = 3000 for different values of the 
transverse size L y of our simulation box, in a flow with A = 640c/w p , a = 10, and 
a — 0.1. We vary L y from L y — 800 c/oj p down to L y — 12.5 c/w p , which approaches 
the ID setup of Petri & Lyubarsky (2007). The black line in the subpanel shows the 
average downstream Lorentz factor (7) as a function of the box width L y . 



the sheet, and its size approaches L y . At this point, a second shock forms, located 
behind the fast shock. For L y > 400 c/u p , this would correspond to the hydrodynamic 
shock discussed in §3, where the striped structure of the flow is erased, and field energy 
is entirely transferred to the particles. Instead, in the case L y < 300c/w p , complete 
dissipation of the alternating fields is prohibited by the fact that the growth of islands 
is artificially inhibited by the transverse extent of the box. Here, the post-shock flow 
retains a striped structure, with regions of hot plasma separated by highly magnetized 
walls of cold particles. With decreasing L y , less magnetic energy is transferred to 
the particles behind the fast shock, and the post-shock flow gets more dominated by 
magnetic field energy, with respect to particle kinetic energy (see the decrease in average 
particle energy with decreasing L y , in the subpanel of Fig. Bl). This also explains why 
in Fig. Bl the normalization of the low-energy peak, populated by the cold particles 
in the high-field regions, grows with decreasing L y , at the expense of the high-energy 
component of hot particles that gained energy from field dissipation. 

As the downstream region gets more magnetized with decreasing L y , the main 
shock, which would be moving at /3 s h ~ 1/3 for L y > 400c/o;p, propagates at a faster 
velocity, eventually catching up with the fast MHD shock, in the limit L„ < A. In this 
regime, we recover the ID results by Petri & Lyubarsky (2007). For the parameters 
employed here, the ID model by Petri & Lyubarsky (2007) would predict negligible 
field dissipation, in agreement with our results for a very narrow box (yellow line for 
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L y = 12.5c/oj p ; see also the subpanel, for small values of L y ). However, as clarified 
by Fig. Bl, this is just an artificial consequence of the reduced dimensionality of their 
model, which cannot correctly capture the development of the tearing-mode instability, 
and the resulting growth and coalescence of magnetic islands. In a more realistic 
2D scenario, complete field dissipation is achieved by the time the flow enters the 
hydrodynamic shock. This stresses that multi-dimensional simulations are essential 
for our understanding of shock-driven reconnection. 
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